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Abstract We study how correlations in the random fitness assignment may aifect the structure of 
fitness landscapes. We consider three classes of fitness models. The first is a continuous phenotype 
space in which individuals are characterized by a large number of continuously varying traits such as 
size, weight, color, or concentrations of gene products which directly affect fitness. The second is a 
simple model that explicitly describes genotype-to-phenotype and phenotype-to-fitness maps allowing 
for neutrality at both phenotype and fitness levels and resulting in a fitness landscape with tunable cor- 
relation length. The third is a class of models in which particular combinations of alleles or values of 
phenotypic characters are "incompatible" in the sense that the resulting genotypes or phenotypes have 
reduced (or zero) fitness. This class of models can be viewed as a generalization of the canonical Bateson- 
Dobzhansky-Muller model of speciation. We also demonstrate that the discrete NK model shares some 
signature properties of models with high correlations. Throughout the paper, our focus is on the perco- 
lation threshold, on the number, size and structure of connected clusters, and on the number of viable 
genotypes. 
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1 Introduction 



The n otion of fitness l andscapes, introduced by a theoretical evolutionary biologist Sewall Wright in 
19321 (see also iKauffman 1993 : Gavrilets .2004, ). has proved extremely useful both in biology and well 
outside of it. In the standard interpretation, a fitness landscape is a relationship between a set of genes (or 
a set of quantitative characters) and a measure of fitness (e.g. viability, fertility, or mating success). In 
Wright's original formulation the set of genes (or quantitative characters) is the property of an individual. 
However, the notion of fi tness landscapes can be generalized to the level of a mating pair, or even a 
population of individuals (|Gavriletsl . l2004h . 

To date, most emp irical information on fitness landscapes in biological applications has come from 
studies of RNA (e.g.. ISchusterl 119951: iHuvnen et all 199a : Fontana and Schusten 119981). proteins (e.g ., 
Lipman and W ilbur"l99l': 'Martinez et al." 1996': 'Rost" 1997'). viruses (e.g..lBurch and Chao"l999'. 'lOm, 
bacteria (e.g.. lElena and Lenski ,2003, : Woods et al, ,2006.) . and artificial life (e.g., Lenski et al.. ,19991 : 
WiUce et al.ll200lh . The three paradigmatic landscapes — rugged, single-peak, and fiat — emphasizing 
particular features of fitnes s landscapes hav e been the focus of most of the earlier theoretical work (re- 
viewed in lKauffmanlll993l : lOavriletsI 120041 ). These landscapes have found numerous app lications with 
regards to the dynamics of adaptatio n (e.g., [K auff man and Levinlll987l : lKauffmanlll993l : iQrri 



and neutral molecular evolution (e.g.. lDerrida and Peliti 



1991 ). 



More recently, it was realized that the dimensionality of most biologically interesting fitness land- 
scapes is enormous and that this huge dimensionality brings some new properties which one does not 
observe in low-dimensional landscapes (e.g. in two- or three-dimensional geographic landscapes). In 
particular, multidimensional landscapes are generically characterized by the existence of neutral and 
nearly neutral networks (also referred to as holey fitness landscapes) that extend throughout the land- 
scapes and that can dramatically affect the evol i itionary dynamics of the populations dGavriletsl. Il997l: 
Gavrilets and Gravnei] . [l997l : lReidvs et al.lll997l:[Gavriletsl . ]2004l : lReidvs etalll200ll : lReidvs and Stadleil . 
2001,.,2002l\ 
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An important property of fitness landscapes is their correlation pattern. A common measure for the 
strength of dependence is the correlation function p measuring the correlation of fitnesses of pairs of 
individual at a distance (e.g., Hamming) d from each other in the genotype (or phenotype) space: 

cov[w(.),w(.)]rf 

pw = r~\ — 



(|Eigen et al.L 119891) . Here, the term in the numerator is the covariance of fitnesses of two individuals 
conditioned on them being at distance d, and vai{w) is the variance in fitness over the whole fitness 
landscape. For uncorrelated landscapes, p{d) = for d > 0. In contrast, for highly correlated land- 
scapes, p{d) decreases with d very slowly. 



The aim of this paper is to extend our previous work (iGavrilets and Gravneri. 119971) in a number of 
directions paying special attention to the question of how correlations in the random fitness assignment 
may affect the structure of genotype and phenotype spaces. For the resulting random fitness landscapes, 
we shed some light on issues such as the number of viable genotypes, number of connected clusters of 
viable genotypes and their size distribution, existence thresholds, and number of possible fitnesses. 

To this end, we introduce a variety of models, which could be divided into two essentially different 
classes: those with local correlations, and those with global correlations. As we will see, techniques used 
to analyze these models, and answers we obtain, differ significantly. We use a mixture of analytical and 
computational techniques; it is perhaps necessary to point out that these models are very far from trivial, 
and one is quickly led to outstanding open problems in probability theory and computer science. 



We start (in Section 2) by briefly reviewing some results from IGavrilets and Gravnen (|l997n . In Sec 



tion 3 we generalize these results for the case of a continuous phenotype space when individuals are 
characterized by a large number of continuously varying traits such as size, weight, color, or the con- 
centrations of some gene products. The latter interpretation of the phenotype space may be particularly 
relevant given the rise of proteomics and the growing interest in gene regulatory networks. 

The main idea behind our local correlations model studies in Section 4 is fitness assignment con- 
formity. Namely, one randomly divides the genotype space into components which are forced to have 
the same phenotype; then, each different phenotype is independently assigned a random fitness. This 
leads to a simple two-parameter model, in which one parameter determines the density of viable geno- 
types, and the other the correlations between them. We argue that the probability of existence of a giant 
cluster (which swallows a positive proportion of all viable genotypes) is a non-monotone function of the 
correlation parameter and identify the critical surface at which this probability jumps almost from to 
1. In Section 4 we also investigate the effects of interaction between conformity structure and fitness 
assignment. 

Section 5 introduces our basic global correlation model, one in which genotypes are eliminated due 
to random pairwise incompatibilities between alleles. This is equivalent to a random version of SAT 
problem, which is the canonical constraint satisfaction problem in computer science. In general, a SAT 
problem involves a set of Boolean variables and their negations that are strung together with OR sym- 
bols into clauses. The clauses are joined by AND symbols into a formula. A SAT problem asks one 
to decide, whether the variables can be assigned values that will make the formula true. An important 
special case, K-SAT, has the length of each clause fixed at K. Arguably, SAT is the most important class 
of problems in complexity theory. In fact, the general SAT was the first known NP-complete problem 
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and was established as such by S. Cook in 1971 (ICooklll97 ih . Even considerable simplifications, such 
as the 3- SAT (see Section 5.4), remain NP- complete, although 2- S AT (see Section 5.1) can be solved 
efficiently by a simple algorithm. See e.g. ^ Korte and Vygen (l2005 h for a comprehensive presentation 
of the theory. Difficulties in analyzing random SAT problems, in which formulas are chosen at random, 
in many ways mi r ror their complexity classes, but even random 2-SAT presents significant challenges 
dde la Vegall2001uBollobas et al.lll994r) . In our present interpretation, the main reason for these difficul- 
ties is that correlations are so high that the expected number of viable genotypes may be exponentially 
large, while at the same time the probability that even one viable genotype exists is very low. In Section 
5, we further illuminate this issue by showing that connected viable clusters must contain fairly large 
sub-cubes, and that the number of such clusters is, in a proper interpretation, finite. The relevance to 
both types of models for discrete and continuous phenotype spaces is also discussed, with particular em- 
phasis on the existence of viable phenotypes in the presence of incompatibihties. Section 5 also contains 
a brief review of the existing theory on higher order incompatibilities. 

In Section 6 we demonstrate how the discrete NK model shares some signature properties of models 
with high correlations. In Section 7 we summarize our results and discuss their biological relevance. The 
proofs of our major results are relegated to Appendices A-E. 



2 The basic case: binary hypercube and independent binary fitness 



We begin with a brief review of the basic setup, from lGavrilets and Gravner (Il997h and lOavriletsI (l2004l\ 



The binary hypercube consists of all n-long arrays of bits, or alleles, that is ^ = {0, 1}". This is our 
genotype space. Genotypes are linked by edges induced by bit-flips, i.e., mutations at a single locus, for 
example, for n = 4, a sequence of mutations might look like 

0000 ^ 1000 ^ 1001 ^ 1101 ^ 1100. 

The (Hamming) distance d{x, y) between x ^ Q and y € ^ is the number of coordinates in which x and 
y differ or, equivalently, the least number of mutations which connect x and y. 

The fitness of each genotype x is denoted by w{x). We will describe several ways to prescribe the 
fitness w at random, according to some probability measure P on the 2^" possible assignments. Then we 
say that an event An happens asymptotically almost surely (a. a. s.) if P{An) ^ 1 as n ^ oo. Typically, 
An will capture some important property of (random) clusters of genotypes. 

We commonly assume that w(x) g |0, 1| so that x i s eithe r viable (w{x) = 1) or inviable {w{x) = 



0). As a natural starting point, iGavrilets and Gravned (| 19971 ) considered uncorrected landscapes, in 



which w{x) is chosen to be 1 with probability p.^, for each x independently of others. We assume this 
setup for the rest of this section and note that this is a well-studied problem in mathematical literature, 
although it presents considerable technical difficulties and some issues are still not completely resolved. 

Given a particular fitness assignment, viable genotypes form a subset of Q, which is divided into 
connected components or clusters. For example, with n = 4, if 0000 is viable, but its 4 neighbors 1000, 
0100, 0010, and 0001 are not, then it is isolated in its own cluster. 

Perhaps the most basic result determines the connectivity threshold ( Toman . 1979 *): when py > 1/2, 



the set of all viable genotypes is connected a. a. s. By contrast, when p^ < 1/2, the set of viable 
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genotypes is not connected a. a. s. This is easily understood, as the connectedness is closely linked to 
isolated genotypes, whose expected number is 2"pt,(l — Pv)""- This expectation makes a transition from 
exponentially large to exponentially small at py = 1/2. The events {x is isolated}, x € C/, are only 
weakly correlated, which implies that when < 1/2 there are exponentially many isolated genotypes 
with high probability, while when p^ > 1/2, a separate argument shows that the event that the set of 
viable genotypes contains no isolated vertex but is not connected becomes very unlikely for large n. 
This is perhaps the clearest instance of the local method: a local property (no isolated genotypes) is 
a. a. s. equivalent to a global one (connectivity). 

Connectivity is clearly too much to ask for, as p^ above 1/2 is not biologically realistic. Instead, 
one should look for a weaker property which has a chance of occurring at small py. Such a property is 
percolation, a. k. a. existence of the giant component. For this, we scale p^ = Xy/n, for a constant A^. 
When > 1, the set of viable genotypes percolates, that is, it a. a. s. contains a component of at least 
c • n~^2" genotypes, with all other components of at most polynomial (in n) size. When A^, < 1, the 
large st component is a. a. s . of size Cn. Here and below, c and C are some constants. These are results 
from lBollobas et all (Il994) . 

The local method that correctly identifies the percolation threshold is a little more sophisticated than 
the one for the connectivity threshold, and uses branching processes with Poisson offspring distribu- 
tion — hence we introduce notation Poisson(A) for a Poisson distribution with mean A. Viewed from, 
say, genotype ... 0, the binary hypercube locally approximates a tree with uniform degree n. Thus 
viable genotypes approximate a branching process in which every node has the number of successors 
distributed binomially with parameters n — 1 and p, hence this random number has mean about A^, and 
is approximately Poisson(At,). When A^, > 1, such a branching process survives forever with probability 
1 — 6 > 0, where 6 = 6{Xy), and 6{X) is given by the imphcit equation 

5 = e^^'-'\ (2) 

(e.g., Athreya and Neylfl97l b. Large trees of viable genotypes created by the branching processes which 



emanate from viable genotypes merge into a very large ("giant") connected set. On the other hand, when 
A„ < 1 the branching process dies out with probability 1. 

The condition A^ > 1 for the existence of the giant component can be loosely rewritten as 



Pi 



1 

> -. 

n 



(3) 



This shows that the larger the dimensionality n of the genotype space, the small er values of the proba- 
bility o f being viable p,, will re s ult in the existe n ce of the giant c omponent. See iGavrilets and Gravner 
(Il997h : lOavriletsI d 19971 . \2004) : Iskippej t.004) : IPigUuccil (Eooeh for discussions of biological signifi- 
cance and implications of this important result. 



3 Percolation in a continuous phenotype space 

In this section we will assume that individuals are characterized by n continuous traits (such as size, 
weight, color, or concentrations of particular gene products). To be precise, we let V = [0, 1]" be the 
phenotype space. 
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We begin with the extension of the notion of independent viability. The most straightforward ana- 
logue of the discrete genotype space considered in the previous section involves Poisson point location 
in V, obtained by generating a Poisson(A) random variable N, and then choosing points xi, . . . , G V 
uniformly at random. These will be interpreted as peaks of equal height in the fitness landscape. Another 
parameter is a small r > 0, which can be interpreted as measuring how harsh the environment is: any 
phenotype within r of one of the peaks is declared viable and any phenotype not within r of one of the 
peaks is declared inviable. For simplicity, we will assume "within r" to mean that "every coordinate 
differs by at most r," i.e., distance is measured in the (n-dimensional) l°° norm || • ||oo- Note that this 
makes the set of viable genotypes correlated, albeit the range of correlations is limited to 2r. 

Our most basic question is whether a positive proportion of viable phenotypes is connected together 
into a giant cluster. Note that the probability that a random point in V is viable is equal to the 
probability that there is a "peak" within r from this point. Therefore, 

p^ = l- exp [-A(2r)"] « A(2r)". 

This is also the expected combined volume of viable phenotypes. 

We will consider peaks xi and Xj to be neighbors if they share a viable phenotype, that is, if their 
r-neighborhoods overlap, or equivalently, if — Xj||oo < 2r. Two viable phenotypes yi and 1/2 are 
connected if they are, respectively, within r of peaks xi and X2, and xi and X2 are connected to each 
other via a chain of neighboring peaks. 

By the standard branching process comparison, the necessary condition for the existence of a giant 
cluster is that a "peak" x is connected to more than one other "peak" on the average. All peaks within 2r 
of the focal peak are connected to the latter. Therefore, if p is the expected number of peaks connected 
to X, then 

// = A-(4rr, 

and ^ > 1 is necessary for percolation. As demonstrated by [Penrose (Il996 h (for a different choice of 



the norm, but the proof is the same), this condition becomes sufficient when n is large. Note that the 
expected number A of peaks can be written as • (4r)^". 

If /i > 1 and fixed, then a. a. s. a positive proportion of all peaks (that is, cN peaks, where c = 
c{p) > 0) are connected in one "giant" component, while the remaining connected components are all 
of size 0(log N). On the other hand, if ^ < 1, all components are a. a. s. of size 0(log N). 

The condition ^ > 1 for the existence of the giant component of viable phenotypes can be loosely 
rewritten as 

Pv > ^. (4) 

This shows that viable phenotypes are likely to form a large connected cluster even when one is very 
unlikely to hit one of them at random, if n is even moderately large. The same conclusion and the same 
threshold are valid if instead of n-cubes we use n-spheres of a constant radius. 

The percolation threshold in the continuous phenotype space given by inequality ([Hi is much smaller 
than that in the discrete genotype space which is given by inequality An intuitive reason for this is 
that continuous space offers a viable point a much greater opportunity to be connected to a large cluster. 
Indeed, in the discrete genotype space there are n neighbors per each genotype. In contrast, in the 



4 PERCOLATION IN A CORRELATED LANDSCAPE WITH PHENOTYPIC 
NEUTRALITY 7 



continuous phenotype space, the ratio of the volume of the space where neigboring peaks can be located 
(which has radius 2r) to the volume of the focal n-cube (which has radius r) is 2". 



4 Percolation in a correlated landscape with phenotypic neutrality 

The standard paradigm in biology is that the relationship between genotype and fitness is mediated 
by phenotype (i.e., observable characteristics of individuals). Both the genotype-to-phenotype and 
phenotype-to-fitness maps are typically not one-to-one. Here, we formulate a simple model capturing 
these properties which also results in a correlated fitness landscape. Below we will call mutations that do 
not change phenotype conformist. These mutations represent a subset of neutral mutations that do not 
change fitness. 

We propose the following two-step model. To begin the first step, we make each pair of genotypes 
X and y in a binary hypercube Q independently conformist with probability Pd(x,y) where d{x, y) is 
the Hamming distance between x and y. We then declare x and y to belong to the same conformist 
cluster if they are linked by a c hain of conformist pairs. This version of long-range percolation model 
(cf., feerger 2004 : Biski^ 2004) divides the set of genotypes Q into conformist clusters. We postulate 



that all genotypes in the same conformist cluster have the same phenotype. Therefore, genetic changes 
represented by a change from one member of a conformist cluster to another (i.e., single or multiple 
mutations) are phenotypically neutral. 

In the second step, we make each conformist cluster independently viable with probability p„ = 
A„/n. This generates a random set of viable genotypes, and we aim to investigate when this set has a 
large connected component. 

For example, the "genotype" can be a linear RNA sequence. This sequence folds into a 2-dimensional 
molecule which has a particular structure (or "shape"), and corresponds to our "phenotype." Finally, the 
molecule itself has a particular function, e.g., to bind to a specific part of the cell or to another molecule. 
A measure of how well this can be accomplished is represented by our "fitness." 

The distribution of conformist clusters depends on the probabilities pi,p2,P3, ■ ■ ■ which determine 
how the conformity p robabihty varies w ith distance. Here we will study the case when pi = pe > 



0,p2 = P3 = ■■■ = (iHaggstroml 120011) . It is then very convenient for the mathematical analysis that 
a pair x and y can be conformist only when they are linked by an edge — therefore we can talk about 
conformist edges or equivalently conformist mutations. (Note however that it is possible that nearest 
neighbors x and y are in the same conformist cluster even if the edge between them is non-conformist.) 

Figure 1 illustrates our 2-step procedure on a four-dimensional example. 

We expect that a more general model with pi declining fast enough with i is just a smeared version of 
this basic one, and its properties are not likely to differ from those of the simpler model. We conjecture 
that for our purposes, "fast enough" decrease should be exponential with a rate logarithmically increasing 
in the dimension n, e.g. for large k, 

Pk < exp(-a(logn)A;), 

for some a > 1. (This is expected to be so because in this case the expected number of neighbors of the 
focal genotype is finite.) 
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Figure 1: A four-dimensional example: start with the cube (top left), create conformist clusters by 
randomly eliminating each edge with probability 1 —pe (top right), remove each conformist cluster with 
probability 1 — Pv (bottom left, removed vertices are black) and finally consider connected components 
of the remaining vertices (bottom right, there is just one component in this case). 



We observe that the first step of our procedure is an edge versio n of the percolation m odel discussed 



in the second section, with a similar giant component transition (jBollobas et all Il992h . Namely, let 
Pi = Pe = Ae/n. Then, if Ae > 1, there is a. a. s. one giant conformist cluster of size c • 2", with all 
others of size at most Cn. In contrast, if Ae < 1 all conformist clusters are of size at most Cn. Note that 
the number of conformist clusters is always on the order 2". In fact, even the number of "non-conformist" 
(i.e., isolated) clusters is a. a. s. asymptotic to e~'^'=2", as P{x is isolated) = (1 — Ae/n)". 

Denote by x <^ y (resp. x ^ y) the event that x and y are (resp. are not) in the same conformist 
cluster. First, we note that the probability P{x y) that two genotypes belong to the same conformist 
cluster depends on the Hamming distance d{x,y) between them, and on = Ag/n. In particular, we 
show in Appendix A that, if Ae < 1 and d{x, y) = kis fixed, then 



klp'^il - Oin'^)) <P{x^y) < A;!p^(l + 0(n"^ log n)). 



(5) 



The dominant contribution k\p^ is simply the expected number of conformist pathways between x and y 
that aie of shortest possible length. 

It is also important to note that, for every x £ Q, the probability P{x is viable) = p^, therefore it 
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does not depend on pe- Moreover, for x, y G G, 
P{x and y viable) — pf, 

= P{x and y viable, x y) + P{x and y viable, x ^ y) — p^ 
= PvP(x ^y)+pl- P{x *^y) -pI 
= p^{l - Pv)P{x <^ y)>0. 

Therefore, the correlation function ([T|l is 

p{x,y) = P{x ^ y), (6) 

which clearly increases with Pf. and, thus, with Ag. Therefore, this model has tunable positive correlations 
controlled by the parameter Ae, whose value does not affect the expected number of viable genotypes. 
The correlation function p{x,y) decreases exponentially with distance d{x,y) when Ag < 1, and is 
bounded below when Ag > 1. Nevertheless, as we will see below, we can effectively use local methods 
for all values of Ag. 



4.1 Threshold surface for percolation 

Proceeding by the local branching process heuristics, we reason that a surviving node on the branching 
tree can have two types of descendants: those that are connected by conformist mutations and those that 
are in different conformist clusters and thus independently viable. Therefore the number of descendants 
is approximately Poisson(Ae + A„). This can only work when Ag < 1, as otherwise the correlations are 
global. 

If Ae > 1, we need to eliminate the entire conformist giant component, which is a. a. s. inviable. 
Locally, we condition on the (supercritical) branching process of the supposed descendant to die out. 
Such conditioned process is a subcritical branching process, with Poisson {XeS) distribution of successors 
(lAthreya and Neyl Il971r) where 6 = 6{Xe) is given by the equation (O. This gives the conformist 
contribution, to which we add the independent Poisson(A„(5) contribution. 

To have a convenient summary of the conclusions above, assume that Ag is fixed and let C('^e) be the 
smallest A^, which a. a. s. ensures the giant component, i.e., 

C{Xe) = inf{A^, : a cluster of at least cn^^2" viable genotypes exists a. a. s. for some c > 0}. 

One would expect that for A^, < ({Xe) all components are a. a. s. of size at most Cn. The asymptotic 
critical curve is given by A^, = C(^e), where 

ifAe[0,l], 
if A G [l,oo). 




Having only a heuristic proof of this, we resort to computer simulations for confirmation. For this, 
we indicate global connectivity with the event A that a genotype within distance 2 of ... is connected 
(through viable genotypes) to a genotype within distance 2 of 1 ... 1. We make this choice because the 
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Figure 2: Simulated A"^ (long dashes) and A^^ (short dashes), and Q (solid) plotted against Ag, for 
n = 10, . . . , 20, and models from Section 4.1 (left frame) and Section 4.2 (right frame). Lower bounds 
increase with n, and upper bounds decrease, for this range of n. 

distance 2 is the smallest that works with asymptotic certainty. Indeed, the genotypes ... and 1 ... 1 
are likely to be inviable. Even the number of viable genotypes within distance one of each of these is 
only of constant order, so even in the percolation regime the probability of connectivity between a viable 
genotype within distance one of ... and a viable one within distance one of 1 ... 1 does not converge 
to 1 but is of a nontrivial constant order. By contrast, there are about r? vertices within distance 2 of 
. . . among which of order n ai^e viable. 

When A^ > C('^e) the probability of the event A should therefore be (exponentially) close to 1. On 
the other hand, when A„ < C(^e) the probability that a connected component within distance 2 of either 
... or 1 ... 1 extends for distance of the order n is exponentially small. We further define the critical 
curves 

= the smallest A^, for which P[A) > 0.1, 
A^ = the largest A„ for which P(^) < 0.9. 

We approximated A^ and A^ forn = 10, . . . , 20 and Ae = 0(0.1)2, with 1000 indepe ndent real- 



ization s of each choice of n, Ae, and A^. We used the linear cluster algorithm described in Sedgewick 



(119971) . The results are depicted in Figure 1. Unfortunately, simulations above n w 20 are not feasible. 
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From Figure 2 we observe that: 

• Even for low n, both critical curves approximate well the overall shape of the theoretical limit 
curve C- 

• and \^ get closer faster than they converge to (. Consequently, one can expect that P{A) 
makes a very sharp jump from near to near 1 even for moderate n. 

• For Ae < 1, A™ tends to be above the Umit curve. This is not really surprising, as the local argument 
always gives an upper bound on the probability P{A) of event A. Further, the approximation of 

deteriorates near Ae = 2, which stems from the possibility of survival of the giant component 
in this regime. 

What is clear from the heuristics and simulations is that conformist mutations, and thus correlations, 
significantly affect the probability of long range connectivity in the genotype space. The effect is not 
monotone: the most advantageous choice is when the correlations are at the point of phase transition 

between between local and global. 

To intuitively understand why percolation occurs the easiest with Ae ~ 1, it helps to think of the 
model as a branching process on clusters rather than on genotypes. For a genotype on a viable cluster, 
there is a number of neighboring clusters and each of these is viable with probability pv If Ae < 1, 
then the probability that any two of the neighboring genotypes are in the same cluster is o(l), so there 
are asymptotically exactly n clusters neighboring the present cluster. Consequently, the overall number 
of descendants will be greater if the size of these clusters is greater on average; which is exactly what 
happens as Ag increases towards 1. If Ag > 1, then there is a positive proportion of the neighboring 
genotypes that are in the giant cluster. This giant cluster is likely to be inviable, so the parameter A^ must 
be greater to compensate for its loss. 

4.2 Correlations between conformity and viability 

In the previous model, the viability probability was independent of the conformity structure. Mainly 
to investigate the robustness of our conclusions, we consider a simple generalization in which there are 
either positive or negative correlations between conformity and fitness. While more sophisticated models 
are possible, the one below is chosen for its amenability to relatively simple analysis. 

Assume now that conformist clusters are formed as before (i.e., with edges being conformist with 
probability pe = Ae/n), are still independently viable, but now the probability of their viability de- 
pends on their size. We will consider the simple case when an isolated genotype (one might call it 
non-conformist) is viable with probability po = Xq/ n, while a conformist cluster of size larger than 1 is 
viable with probability p\ = Xi/n. 

In this case 

P{x is viable) = (1 -pe)>o + (1 - (1 -PeT)pi - ^ (e-^^Ao + (1 - e-^^)Ai) . 
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Moreover, by a similar calculation as before, 

P{x and y viable) — P{x viable)^ 

= - pi)P{x ^y) + P{x non-conformist )2pe(po - Pif ■ l{d{x,y)=i}- 

Here, the last factor is the indicator of the set {{x, y),d{x, y) = 1}, which equals 1 if d{x, y) = 1 and 
otherwise. Therefore, for d{x, y) > 2, the correlation function ([T]) is 

~ e-A.Ao + a-e'A.);,/(- ^ v)' 
which is smaller than before iff Ai < Aq- However, it has the same asymptotic properties unless Ai = 0. 



Assume first that Ae < 1. The local analysis now leads to a multi-type branching process (lAthreya and NeyL 



Il97lh with three types: NC (non-conformist node), CI (non-isolated node independently viable, so no 
conformist edge is accounted for), and CC (non-isolated node viable by conformity, so a conformist edge 
is accounted for). 

Note first that a genotype is non-conformist with probability about 6^^"=. Hence a node of any of 
the three types creates a Poisson(e^'*''= Ai) number of type NC descendants, and a Poisson((l — e~'^'=)Ai) 
number of type CI descendants. In addition, the type CI creates a Poisson(Ae), conditioned on being 
nonzero, number of descendants of type CC and type CC creates a Poisson(Ae) number of descendants 
of type CC. Thus the matrix of expectations, in which the ijth entry is the expectation of the number of 
type j descendants from type i, is 



M 



e ^=Ao (l-e-^^)Ai Ae/(1 - e"^^ 



= Ao (l-e-^^)Ai Ae 
When Ae > 1, Ae needs to be replaced by Xe6, and Ai by Xi6, where 6 = 6{Xe) is given by 



It follows from the theory of multi-type branching processes (|Athreya and Neyl.ll97Ih that the critical 
surface for survival of a multi-type branching process is given by det(M — 1) = 0. 

The simplest case is when only non-conformist genotypes may be viable, i.e., Ai = 0. In this 
case the critical surface is given by Xoe~'^^ = 1 (Pitman, unpub.). Not surprisingly, the critical Aq to 
achieve global connectivity strictly increases with Ae, which is the result of negative correlations between 
conformity and viability. 

The other extreme is when non-conformist genotypes are inviable, i.e., Aq = 0. As an easy compu- 
tation demonstrates, the critical curve is now given by Ai = Ci^e), where 

C(A) = I Ae-'+l-e- if AG {0,1}, 

Note that C(A) — > oo as A ^ 0. We carried out exactly the same simulations as before. These are 
also featured in Figure 2 (right frame), and again confirm our local heuristics. We conclude that positive 
correlations between viability and conformity tend to lead to a V-shaped critical curve, whose sharpness 
at critical conformity Ae = 1 increases with the size of correlations. In short, then, correlations help 
more if viability probability increases with size of conformist clusters. 
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5 Percolation in incompatibility models 



In the model considered in the previous section correlations rapidly decreased with distance. This prop- 
erty made local analysis possible. The models we introduce now are fundamentally different in the sense 
that correlations are so high that the local method gives a wrong answer. 

In the previous sections, in constructing fitness landscapes we were assigning fitness to individ- 
ual genotypes or phenotypes. Here, we make certain assumptions about "fitness" of particular com- 
binations of alleles or the values of phenotypic characters. Specifically, we will assume that some of 
these combinations are "incompatible" in the sense that the resul t ing ge notypes or phenotypes have re- 
duced (or zero) fitness (On, 1995 : Orr and On . 19961 : Gavrilets . 2004 ) . The resii l ting models can be 



viewed as a generalization of the Bateson-Dobzhansky-MuUer model (lOrri 1 19951: lOrr and Orn. 11996 : 
Orr', '1997'; Orr and Turelli,'200l'; Gavrilets and Hastings. 1996 : Gavrilets . 19971 : Gavrilets and Gravnei . 



1997. : ,Gavrilet&. ,2003i. ,2004 : .Coyne and Orr. .2004.) which represents a canonical model of speciation. 



5.1 Diallelic loci 

We begin by assuming that viability of a genotype is determined by a set F of pairwise incompatibilities. 
F is thus a subset of 4 • (") pairs {ui,Vj), where I < i < j < n and n, v € {0, 1}. In this nonstandard 
notation, (61,02) S F, for example, means that allele at locus 1 and allele at locus 2 are incompatible. 
In general, if {ui ,Vj) € F, all genotypes with u in position i and v in position j are inviable. A genotype 
X is then inviable if and only if there exist i and j, with i < j, so that u and v are, respectively, the 
alleles of x at loci i and j, and (uj, vj) G F. For example, if Fi = {(Oi, O2), (I21 03), (li, I2)}, viable 
genotypes may have Oil, 100, and 101 as their first three alleles. For F2 = Fi U {(Oi, I3), (li, O2)}, no 
viable genotype remains. 

Incompatibility (Oi , O2) is equivalent to two implications: Oi =^ I2 and O2 =^ li or to the single 
OR statement li OR I2. In this interpretation, the problem of w hether, for a given l ist of incompatibilities 



F, there is a viable genotype is known as the 2-SAT problem (IKorte an d Vygenl. 120051) . The associated 



digraph Dp is & graph on 2n vertice s x,:, i = 1, . . . ra, x = , 1, with oriented edges determined by the 



implications. A well-known theorem (IKorte and Vygenll2005h states that a viable genotype exists iff Dp 



contains no oriented cycle from Oj to Ij and back to 0^ for any i = 1, ... n in Dp. For example, for the 
incompatibilities F2 as above, one such cycle is Oi — > I2 — I3 — > li ^ I2 ^ Oi. 

Now assume that each possible incompatibility is adjoined to F at random, independently with prob- 
ability 

c 

^ ~ 2^' 

(We use the generic notation p for a probability parameter in all our models, even though the nature of 
probabilistic assignments differs from model to model.) 

Existence of viable genotypes. Let be the number of viable genotypes. Then 
• if c > 1, then a. a. s. A = 0. 



• if c < 1, then a. a. s. A > 0. 
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This result first appeared in the computer science literature in the 90's (see Ide la Vega |200 ll for a review) 



2001: Janson et al- 



and it is an extension of the celebrated Erdos-Renyi random graph results (iBoUobasl . 
20001) to the oriented case. 

Note that the expectation E{N) = 2"(1 - p)^^) « 2"e-'="/^ which grows exponentially whenever 
c < 4 log 2 2.77. Negle cting correlation s would therefore suggest a wrong threshold for > 0. The 
local method (e.g., used in iGavriletsI l2004l . Chapter 6) is even farther off, as it suggests an a. a. s. giant 
component when p < (1 — e) log n/n for any e > 0. 

The number of viable genotypes. Ass ume that c < 1. Sophis t icated , but not mathema tically 
rigorous methods based on replica symmetry ( Monasson and Zecchina . 1997 : Biroli et al. , 2000h from 
statistical physics suggest that, as n — > oo, limn~^ log A'' varies almost linearly between log 2 ^ 0.69 
(for small c, when, as we prove below, this limit is log 2 + 0{c)) and about 0.38 (for c close to 1). On e 
can however prove that log N is for large n sharply concentrated around its mean (Ide la Vegall200lh . 

Upper and lower bounds on N can also be obtained rigorously. For example, if X is a number of in- 
compatibilities which involve disjoint pairs of loci (i.e., those for which every locus is represented at most 
once among the incompatibilities), then N < exp(n log 2 + X log(3/4)), as each of the X incompatibil- 
ities reduces the number of viable genotypes by the factor 3/4. If we imagine adding incompatibilities 
one by one at random until there are about cn of them, then after we have k incompatibilities on disjoint 
pairs of loci the waiting time (measured by the number of incompatibilities added) for a new disjoint 



one is geometric with expectation {^) / 
approximate equation 



Therefore, X is a. a. s. at least Kn, where K solves the 




~ cn, 



or 



Kn 



1 



dk ^ —, 

n 



which reduces to K 



/o (n - 2kf 

c/(l + 2c). This implies that the upper bound on N can be defined as 



lim sup — log N < 

n 1 + 2c 



log 2 + 



1 + 2c 



log 3. 



(9) 



A lower bound is even easier to obtain. Namely, the probability that a fixed location (i.e., locus) i 



.p)4(n- 

- e" 



1) 



e and then it is easy to see that the number of loci represented 



does not appear in F is ( 1 - 

in F is asymptotically (1 — e~'^'^)n. As the other loci are neutral (in the sense that changing their alleles 
does not affect fitness), n^^ log N is asymptotically at least e^'^'^ log 2. Clearly, this gives a lower bound 
on the exponential size of any cluster of viable genotypes. 

If this was an accurate bound, it would imply that the space of genotypes is rather simple, in that 
almost all its entropy would come from neutral loci. The Appendix B presents two arguments which will 
demonstrate that this is not the case. The derivations there are somewhat technical, but do provide more 
insight into random pair incompatibilities. 

The structure of clusters. The derivations in Appendix B show that every viable genotype is 
connected through mutation to a fairly substantial viable sub-cube. In this sub-cube, alleles on at most 
a proportion ru(c) < 1 of loci are fixed (to or 1) while the remaining proportion 1 — ru{c) could be 
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Figure 3: Simulated number of clusters, vs. c for n = 20. The proportion (out of 1000) of trials with 
exactly one, exactly two, and at least three clusters is plotted respectively with +'s, x's and *'s. The 
solid curve is 1^(1 — c)e'^. 



varied without effect on fitness. Note from Figure 4 in the Appendix B that 1 — r„(c) > 0.3 for all c, 
and that such a phenomenon is extremely unlikely on uncorrected landscapes. Note also that, for c < 1, 
N > 2(i^'^'"('=))" a. a. s. and so the lower bound on N can be written as 



liminf — log > (1 — J^ufc)) log 2. 
n 



(10) 



The number of clusters. The natural next question concerns the number of clusters R when c < 1. 
This again has quite a surprising answer, unparalleled in landscapes with rapidly decaying correlations. 
Namely, R is stochastically bounded, that is, for every e > there exists an z = z(e) such that P{R < z 
for all n) > 1 — e. A s there is some co nfusion in the literature as to whether it is even possible to get 
more than one cluster ( Biroli et al. , 2000l) . Appendix C presents a sketch of the results which will appear 



in Pitman (unpub.). There we also show that the limiting probability of a unique cluster is i/(l — c)e'^. 

Asymptotically, a unique cluster has a better than even chance of occurring for c below about 0.9, and 
is very likely to occur for small c, though of course not a. a. s. so. To confirm, we have done simulations 
for n = 20 and c = 0.01(0.01)1 (again 1000 trials in each case) and got distribution of clusters depicted 
in Figure 3. The results suggest that the convergence to limiting distribution is rather slow for c close to 
1 , and that the likelihood of a unique cluster increases for low n. 

To summarize, in the presence of random pairwise incompatibilities, the set of viable genotypes is, 
when nonempty, divided into a stochastically bounded number of connected clusters, where a unique 
cluster is usually the most likely possibility. These clusters are all of exponentially large size (with 
bounds given by equations l9l and [TOl) . in fact they all contain sub-cubes of dimension at least (1 — r„(c))n. 
However, the proportion of viable genotypes among all 2" genotypes is exponentially small, by equation 
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5.2 Multiallelic loci 



Here we assume that at each locus there can be a (> 2) alleles (cf.. iReidysl 120060 . In this case, the 
genotype space is the generalized hypercube = {0, . . . , a — 1}". For a = 3 this could be interp reted as 
the genotype space of diploid organisms without cis-trans effects (|Gavrilets and GravnerL 119971) . o = 4 
corresponds to DNA sequences, and a = 20 corresponds to proteins. Much larger values of a can 
correspond to a number of alleles at a protein coding locus and we will see later that there is not much 
difference between this model and a natural continuous space model. 

We will assume that each pair of alleles, out of total number of a? (2) is independently incompatible 
with probability 

c 

The main question we are interested in here is for which values of c viable genotypes exist a. a. s. 
Clearly, if N is the number of viable phenotypes, then the expectation 

E{N) = a"(l -p)(2) exp(nloga - \cn), 

and so there are a. a. s. no viable phenotypes when c > 4 log a. On the other hand, clearly there are 
viable genotypes (with all positions filled by O's and I's) when c < 1. It turns out that the first of these 
trivial bounds is much closer to the c r itical value when a is large. Before we proceed, however, we state 
a sharp threshold result from lMoUoyl (120031) : there exists a function 7 = 7(77,, a) so that for every e > 0, 



• if c > 7 + e, then a. a. s. = 0. 

• if c < 7 — e, then a. a. s. iV > 0. 

In words, for a fixed a, the probability of the event that > 1 transitions sharply from large to small as 
np varies. As it is not proved that lim„^oo lin, a) exists, it is in principle possible that the place of this 
sharp transition fluctuates as n increases (although it must of course remain within [1, 4 log a]). 

Our main result in this section is 



7 = 4 log a — 0(1), as a ^ 00. 



(11) 



T his somewhat surprising result in p roven in Appendix D by the second moment method, as developed 
in lAchlioptas and Moord (|2004) and lAchhoptas and PeresI (|20041) . 



5.3 Continuous phenotype spaces 

Here we extend the model of pair incompatibilities for the case of continuous phenotypic space V. Again, 
we have a small r > as a parameter. For each of (i, j), i < j, we consider independent Poisson point 
location liij in the unit square [0,1] x [0,1], of rate A = c/(2n). (Equivalently, choose Poisson(A) 
number of points uniformly at random in [0, 1] x [0, 1].) Then we declare a G V inviable if there exist 
i < j so that {ai,aj) is within r of Ylij. Again, we use the two-dimensional norm for distance. 
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Our procedure can be visualized as throwing a random number of (n — 2) -dimensional square tubes of 
inviable phenotypes into the phenotype space. 

Our main result here is that the existence threshold is on the order c f« — log r/r^. Namely, we prove 
in the Appendix E that there exists a constant C > so that for small enough r, 

• if c> 4^^, then a. a. s. iV = 0. 

• if c < then a. a. s. iV > 0. 



5.4 Complex incompatibilities 



Here we assume that incompatibilities involve K {> 2) diallelic loci ( Orr and On . 19961 : Gavrilets . 
2004) . The quest ion whether a viable combination of genes exist is then equivale nt to the K-SAT 
probl em ((Korte and Vygen . 2005 b. Even for K = 3, this is an NP-complete problem ( Korte and VygenL 
20051), so there is no known polynomial algorithm to answer this question. The random case, which we 
now describe, is also much harder to analyze than the 2- SAT one. Let F be a random set to which any 
]^) incompatibilities belong independently with probability 



of the 2^^" 



P 



n 



K-i • 



Here c = c{K) is a constant, and the above form has been chosen to make the number of incompati- 
bilities in F asymptotically cn. (Note also th e agreernent wi th the definition of p in Section 5.1 when 
K = 2.) For a fixed K, it has been proved (|FriedgutL 1 19991 ) that the probability that viable genotype 
exists jumps sharply from to 1 as c varies. However, the location of the jump has not been proved to 
converge as n — » oo. Instea d, a lot of effort has been invested in obtaining good bounds. For example 
( Achlioptas and Peres , 2004 1). for K = 3, c < 3.42 implies a. a. s. existence of viable genotype, while 
c > 4.51 implies a. a. s. nonexistence (while the sharp constant is estimated to be about 4.48, see e.g. 
BiroU et allEoOol) . For K = 4 t he best current bounds are 7 .91 and 10.23. For large K, the transition 
occurs at c = 2^ log 2 - 0{K) ((Achlioptas and PeresLliooi) . 



Techniques from statistical physics (|Biroli et all l2000h strongly suggest that, for K > 3, there is an- 
other phase transition, which for = 3 occurs at about c = 3.96. For smaller c, the viable genotypes are 
conjectured to be contained in a single cluster. For larger c, the space of viable genotypes (if nonempty) 
is divided into exponentially many connected clusters. 

Perhaps more relevant to genetic incornpatibi lities is the following mixed model (commonly known 
as (2 + p)-SAT). iMonasson and Zecchinalll997n . Assume that every 2-incompatibility is present with 
probability C2/(2n), while every 3-incompatibility is present with probability 3c3/(4n^). The normal- 
izations are chosen so that the numbers of the two types of incompatibilities are asymptotically C2n and 
csn, respectively. 

If C2 (resp. C3) is very small, then the respective incompatibility set affects a very small propor- 
tion of loci, therefore C3 (resp. C2) determines whether a viable genotype is likely to exist. Intu- 
itively, one also expects that 2-incompatibilities should be more important than 3-incompatibilities as 
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one of the former type excludes more genotypes than one of the latter type. A careful analysis con- 
firms this. First observe t hat co > 1 implies a. a . s. no n-existence of a viable genotype. The surprise 
(IMonasson and Zecchinal . 119971 : lAchlioptas et all l2001r) is that if C3 is small enough, C2 < 1 implies 
a. a. s. ex istence of viable genotypes, so t he 3 -incompatibilities do not change the thresh old. This is estab 



lished in lMonasson and Zecchinal (119971 ) by a physics argument for C3 < 0.703, while lAchlioptas et al 



(I2OOII) gives a rigorous argument for C3 < 2/3. Therefore, even if their numbers are on the same scale, 
if the more complex incompatibilities are rare enough compared to the pairwise ones, their contribution 
to the structure of the space of viable genotypes is not essential. 



6 Notes on neutral clusters in the discrete NK model 



T he model considered here i s a sp ecial case of the discretized NK model (iKauffmanl . Il993h . introduced 
in 



Newman and Engelhardl (Il998h . This model features n diallelic loci each of which interacts with K 
other loci. To have a concrete example, assume that the loci are arranged on a circle, so that n + 1 = 1, 
n + 1 = 2, etc., and let the interaction neighborhood of the i'th locus consist of itself and K loci to its 
right i + . . . ,i + K. For a given genotype x £ Q = {0, 1}", the neighborhood configuration of the f'th 
locus is then given by Mi{x) = (xj, Xj+i, . . . , Xi+x) G {0, 1}^+^. To each locus and to each possible 
configuration in its neighborhood we independently assign a binary fitness contibution. To be more 
precise, we choose the 2^+^n numbers Vi{y), i = I, . . . ,n and y € {0, 1}^^+^, to be independently or 
1 with equal probability, and interpret Vi{y) as the fitness contribution of locus i when its neighborhood 
configuration is y. The fitness of a genotype x is then the sum of contributions from each locus: 



i=l 



In IKauffmanl (Il993l) . the values Vi were taken from a continuous distribution. In lNewman and Engelhardl 
(| 19981) . these values were integers in the range [0, F — 1] so that our model is a special case F = 2. 
Neutral clusters are connected components of same fitness. 

The K = case is easy but nevertheless illustrative. Namely, a mutation at locus i will not change 
fitness iff Vi{0) = Vi{l); let D be the number of such loci. Then D ~ n/2 a. a. s., the number of different 
fitnesses isn — D, each neutral cluster is a sub-cube of dimension D, and there are exactly 2"^^ neutral 
clusters. 

The next simplest situation is when K = I. Let Di be the number of loci i for which Vi is constant. 
Then Di ~ n/8 a. a. s., and each neutral cluster contains a sub-cube of dimension Di. Moreover, let D2 
be the number of loci i for which (00) = Vi{Ol) / 't^j(lO) = ui(ll). Note that any genotypes that differ 
at such locus i must belong to a different neutral cluster, and so the number of different neutral clusters 
is at least 2^^ xhus there are exponentially many of them, as again D2 ~ n/8 a. a. s. This division of 
genotype space into exponentially many clusters of exponential size persists f or every K, although the 
distri bution of numbers and sizes of these clusters is not well understood (see iNewman and Engelhardl 
19981 for simulations for n = 20). 



Finally, we mention that the question of whether a genotype with the max imal possible fitn ess n 
exists for a given K is in many way related to issues in incompatibilities models dChoi et all 120050 . 
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7 Discussion 

In this section we summarize our major findings and provide their biological interpretation. 

The previous work on neutral and nearly neutral networks in multidimensional fitness landscapes 
has concentrated exclusively on genotype spaces in which each individual (or a group of individuals) is 
characterized by a discrete set of genes. However many features of biological organisms that are actually 
observable and/or measurable are described by continuously varying variables such as size, weight, color, 
or concentration. A question of particular biological interest is whether (nearly) neutral networks are as 
prominent in a continuous phenotype space as they are in the discrete genotype space. Our results provide 
an affirmative answer to this question. Specifically, we have shown that in a simple model of random 
fitness assignment, viable phenotypes are likely to form a large connected cluster even if their overall 
frequency is very low provided the dimensionality of the phenotype space, n, is sufficiently large. In fact, 
the percolation threshold for the probability of being viable scales with n as 1/2" and, thus, decreases 
much faster than 1/n which is characteristic of the analogous discrete genotype space model. 

EarUer work on nearly neutral networks has been limited to consideration of the relationship between 
genotype and fitness. Any phenotypic properties that usually mediate this relationship in real biological 
organisms have been neglected. In Section 4, we proposed a novel model in which phenotype is intro- 
duced explicitly. In our model, the relationships both between genotype and phenotype and between 
phenotype and fitness are of many-to-one type, so that neutraUty is present at both the phenotype and 
fitness levels. Moreover, this model results in a correlated fitness landscape in which the correlation 
function can be found explicitly. We studied the effects of phenotypic neutrality and correlation between 
fitnesses on the percolation threshold and showed that the most conducive conditions for the formation 
of the giant component is when the correlations are at the point of phase transition between local and 
global. To explore the robustness of our conclusions, we then look at a simplistic but mathematically 
illuminating model in which there is a correlation between conformity (i.e., phenotypic neutraUty) and 
fitness. The model has supported our conclusions. 

Section 5, we studied a number of models that have been recentiy proposed and explored within the 
context of studying speciation. In these models, fitness is assigned to particular gene/trait combinations 
and the fitness of the whole organisms depends on the presence or absence of incompatible combinations 
of genes or traits. In these models, the correlations of fitnesses are so high that local methods lead to 
wrong conclusions. First, we established the connection between these models and K-SA1 problems, 
prominent in computer science. Then we analyzed the conditions for the existence of viable genotypes, 
their number, as well as the structure and the number of clusters of viable genotypes. These questions 
have not been studied previously. Among other things we showed that the number of clusters is stochasti- 
cally bounded and each cluster contains a very large sub-cube. The majority of our results are for the case 
of pairwise incompatibilities between diallelic loci, but we also looked at multiple alleles and complex 
incompatibiUties. Moreover, we generalized some of our results to continuous phenotype spaces. 

At the end, we provided some additional results on the size, number and structure of neutral clusters 
in the discrete NK model. 

Some more general lessons of our work are that 

• Correlations may help or hinder connectivity in fitness landscapes. Even when correlations are 



REFERENCES 



20 



positive and tunable by a single parameter, it may be advantageous (for higher connectivity) to 
increase them only to a limited extent. 

• Averages (i.e., expected values) can easily lead to wrong conclusions, especially when correlations 
are strong. Nevertheless, they may still be useful with a crafty choice of relevant statistics. 

• Very high correlations may fundamentally change the structure of connected clusters. For example, 
clusters may look locally more like cubes than trees and their number may be reduced dramatically. 

• Necessary analytical techniques may be unexpected and quite sophisticated; for example, they may 
require detailed understanding of random graphs, spin-glass machinery, or decision algorithms. 
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Appendix 

Appendix A. Proof of equation 

To prove equation (5), we assume that Ag < 1 and show that for a fixed k (which does not grow with 
n), the event that x and y at distance k are in the same conformist cluster is most likely to occur because 
X and y are connected via the shortest possible path. Indeed, the dominant term k\p^ is the expected 
number of conformist pathways between x and y that are of shortest possible length k. This easily 
follows from the observation that on a shortest path there is no opportunity to backtrack; each mutation 
must be toward the other genotype. We can assume that x is the all O's genotype and y is the genotype 
with I's in the first k positions and O's elsewhere. There are kl orders in which the I's can be added. 

To obtain the lower bound we use inclusion-exclusion on the probability that x <^ y through a 
shortest path. Let Ii = Ti{x, y) be the set of all paths of length / between x and y. Then 

P{x ^y)>Y^ P{A^) - ^(^" n "^P) 

where is the event that a particular path a consists entirely of conformist edges. Notice that two 
distinct paths of the same length differ by at least two edges. Thus, we get the following upper bound 



k+2 

'e ) 



and the lower bound in (5) follows. 

The upper bound is a little more difficult to obtain (it is only here that we use Ag < 1) and we 
need some notation. Each genotype can be identified with the set of I's that it contains, so for any two 
genotypes u and u we let it A f denote the set of loci on which they differ. Notice that if n A v is even 
(resp. odd) then every path between u and v is of even (resp. odd) length because each mutation which 
alters the allele at a locus not in u A w must later be compensated for. 

To estimate the expected number of conformist pathways, we will need to bound the number of paths 
of length I between x and y. This is given by 

k\( jmln™" where m = ^. 

\m 2 



We show this via the methods of iBoUobas et al.l (119921) . They obtain an estimate for the number of cycles 



of a given length through a fixed vertex of the cube. 

Given a path, say x = vq.vi, . . . ,vi = y, between x and y, let us associate the sequence (ei^i, . . . , e/i;) 



where 



A. and ^.={!| !f;;^=;;':;\{-^} 



j = Since distinct paths will have distinct sequences we can bound the number of paths by 

finding an upper bound for the number of sequences. 

Note that there must be m + A; positive entries, which occur at (^^^j) = (^) possible locations. The 
absolute values of m of these entries are chosen freely from {1, . . . , n}, while the remaining k must be 
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the integers 1, . . . ,k. There are 77,"*^! ways to do this. We are free to order the m negative entries and 
the bound follows. 

We now assume that d{x, y) is even and relabel d{x, y) = 2k. We omit the similar calculation for 
odd distances. Define b = — 3/c/(21og Ag) and t = [61ognJ. Then the expected number of conformist 
paths between x and y can be expressed as 

E E^" = E E^e' + EE^e' 

i>k+i k+i<i<t i>t 



k+l<l<t ^ ^ l>t 

k+l<l<t l>t 

< {2k)\pf Y {2h\ePe\ogn)^~'^ + 0{\f^°^'') 

l>k+l 

= k{2k)\pfO{pe\ogn) + 0{n^^^°^^^) 
= k{2k)\pfO{n-Hogn). 



Appendix B. Cluster structure under random pair incompatibilities. 



Here we show that, under random pairwise incompatibilities model introduced in Section 5.1, connected 
clusters include large subcubes. The basic idea comes from lBoufkhad and DuboisI (Il999h . A configura- 
tion a G {0, 1, *}" is a way to specify a sub-cube of Q, if *'s are thought of as places which could be 
filled by either a or a 1. The number of non-*'s is the length of a. Call a an implicant if the entire 
sub-cube specified by a is viable. 

We present two arguments, beginning with the one which works better for small c. Let the auxiliary 
random variable X be the number of pairs of loci («, j), i < j, for which: 



(El) There is exactly one incompatibility involving alleles on i and j. 

(E2) There is no incompatibility involving an allele on either i or j, and an allele on A; ^ 

Assume, without loss of generality, that the incompatibility which satisfies (El) is (Ij, Ij). Then fitness 
of all genotypes which have any of the allele assignments OjOj, Ojlj and Ifij, and agree on other loci, is 
the same. Note also that all pairs of loci which satisfy (El) and (E2) must be disjoint. Therefore, if x is 
any viable genotype, its cluster contains an implicant with the number of *'s at least X plus the number 
of free loci. To determine the size of X, note that the expectation 

and furthermore, by an equally easy computation, 



E{X^)-E{Xf = 0{n), 
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so that X ~ ce ^^n a. a. s. It follows that every cluster contains a. a. s. at least exp((e ^^+ce log 2 — 
e)n), viable genotypes, for any e > 0. 

The second argument is a refinement of the one in Boufkhad and DuboisI (1999) and only works 



better for larger c. Call an implicant a a prime implicant (PI) if at any locus i, replacement of either Oj or 
Ij by *i results in a non-implicant. Moreover, we call a the least prime implicant (LPI) if it is a PI, and 
the following two conditions are satisfied. First, if all the *'s are changed to O's, then no change from Ij 
to Oj results in a viable genotype. Second, no change to where i < j, results in an indicator. 

Now, every viable genotype must have an LPI in its cluster. To see this, assume we have a PI for 
which the first condition is not satisfied. Make the indicated change, then replace some O's and I's by *'s 
until you get a prime indicator. If the second condition is violated, make the resulting switch, then again 
make some replacement by *'s until you arrive at a PI. Either of these two operations moves within the 
same cluster, and keeps the number of 1 's nonincreasing and their positions more to the left. Therefore, 
the procedure must at some point end, resulting in an LPI in the same cluster. 

For a sub-cube a to be an LPI, the following conditions need to be satisfied: 

(11) Every non-* has to be compatible with every other non-*, and with both and 1 on each of the *'s. 

(12) Any of the four 0,1 combinations on any pair of *'s must be compatible. 

(LPIl) Pick an i with allele 1, that is, a Ij. Then Oj must be incompatible with at least one non-*, or at 
least one on a *. Furthermore, if Oj has an incompatibility with a on a * to its left, it has to have 
another incompatibility, either with a non-*, or with a or a 1 on a *. 

(LPI2) Pick a Oj. Then Ij must be incompatible with a non-*, or a or a 1 on a *. 

The first two conditions make a an implicant, and the last two an LPI. Note also that these conditions are 
independent. 

Let now X be the number of LPI of length rn. We will identify a function L4 = L4 (r, c) such that 

-logE{X) < L4. 
n 

Let 

Li = Li{P,p, z) = z{(3 logp + (1 - /3) log(l -p)-piogP-{l-P) log(l - 

This is the exponential rate for the probability that in zn Bernoulli trials with success probability p there 
are exactly /?n successes, i.e., this probability is exp(Lin). Further, if K,e,6 G (0, 1) are fixed, then 
among sub-cubes with rn non-*'s and an I's (a < r), the proportion which have en I's in [nn, n] and 
5n *'s in [1, Kn] has exponential rate 

L2 =L2{r,c,K,a,e,5) 

=Li((a — e)/K, a, k) + Li(e/(1 — k), a, 1 — k) 

+ Li((5/(k - a + e), 1 - r, k - a + e) + Li{{l - r-6)/{l - k - e),l - r,l - k - e). 



(Here all four first arguments in Li are in [0, 1], or else the rate is —00.) 
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Figure 4: The upper bound r„(c) for the number of non-*'s in the implicant of smallest length included 
in every cluster of viable genotypes, plotted against c. 

The expected number of LPI, with r, k, e, S given as above, has exponential rate at most (and this is 
only an upper bound) 

L3 =L3{r,c,K,a,€,5) 

= — (1 — r) log(l — r) — a log a — {r — a) Iog(r — a) 
-c{l-r/2f 

+ {r -a) log(l - cxp(-c(l - r/2))) 

+ (a - e) log(l - cxp(-c/2)) + elog(l - exp(-c/2) - i5cexp(-c(l - r/2))) 
+ L2(r, c, K, a, e,5). 

The next to last Une is obtained from (LPIl), as en I's must have 5n *'s on their left. 
It follows that L4 can be obtained by 



If Li{r, c) < 0, all LPI (for this c) a. a. s. have length at most r. Numerical computations show that this 
gives a better bound than 1 — e"^*^ — ce"^*^ for c > 0.38. Let us denote the best upper bound from the 
two estimates by r„(c). This function is computed numerically and plotted in Figure 3. 

Appendix C. Number of clusters under random pair incompatibilities 

In this section we briefly explain why the number of clusters under random pair incompatibilities is 
asymptotically a function of a Poisson random variable. There is a clear way to separate the genotype 
space into disconnected clusters. For example, if Fi = {(Oi, O2), (I2, 03), (li, I2)}, we see that every 
viable genotype has one of these two allele configurations on the first two loci: C = O1I2 or C = I1O2. 
Since there are no genotypes with OiOi or I1I2, there is no way to mutate from the viable genotypes 
with O1I2 to the viable genotypes with I1O2 without passing through an inviable genotype. However, if 



L4(r, c) = inf sup L^{r, c, k, a, e, 6). 

a,e,5 
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we add one incompatibility to Fi to make F2 = FiU{(0i,l2)}, then there are no longer any genotypes 
with the alleles O1I2 and we return to a single cluster of viable genotypes. 

Notice that the digraph Dp-^ contains the directed cycle li ^ O2 — > li and equivalently the directed 
cycle I2 ^ Oi ^ I2. Dp^ also contains these cycles but there are paths between them as well: O2 —>■ Oi 
and li I2. 

Formally, a pair of complementary allele configurations (C, C) on a set of A; > 2 loci is defined to 
be a splitting pair if the digraph Dp contains a directed cycle (in any order) on the alleles in C (and 
equivalently on those in C, which consist of reversed alleles in C) and does not contain a path between 
the alleles in C and the alleles in C. It should be clear from the example Fi above that the existence 
of a splitting pair will create a barrier in the genotype space through which it is not possible to pass by 
mutations on viable genotypes. In fact, it is proved in Pitman (unpub.) that any two viable genotypes u 
and V will be disconnected in the fitness landscape if and only if the loci on which they differ contain a 
splitting pair. 

Thus, the existence of viable genotypes on either side of a splitting pair (with each configuration 
of complementary alleles) ensures disconnected clusters. If there are k splitting pairs in the formula F 
and there are viable genotypes with each of the allele configurations in each of the splitting pairs then 
there are 2^^ clusters of viable genotypes. The restriction that there be viable genotypes on either side 
is asymptotically unlikely to make a difference as we can fix one of the 2^^ configurations of alleles and 
a. a. s. find a viable genotype on the remaining loci. Therefore the number of clusters of viable genotypes 
is a. a. s. equal to 2^, where X is the number of splitting pairs, provided that X is stochastically bounded, 
but we will see shortly that the expectation E{X) is bounded. In fact, the next paragraph suggests that 
X converges to a Poisson limiting distribution. (A detailed discussion of this issue will appear in Pitman 
(unpub.).) 

It follows from IPalastil (Il97lb or lBoUobasI (|200lh that the number of directed cycles of length k in 



Dpi?, Poisson(Afc) with = {2k) ^c^. In particular, the expected number of splitting pairs converges 
to is A = — i(ln(l — c) + c). Moreover, the probability that the re is no splitti ng pair converges to the 



product of the probabilities that the cycle of each length is absent (IPalastil . 1 197 11) . which is 



n 



exp|-|:)=exp(l!i<l^l=[(l-.Kli (12, 



fc=2 

In particular, this gives the limiting probability of a unique cluster. 



Appendix D. Proof of equation (fll]) . 

In this section we assume that genotypes have multiallelic loci, which are subject to random pair incom- 
patibilities. The model introduced in Section 5.2 is the most natural, but is not best suited for our second 
moment approach. Instead, we will work with the equivalent modified model with m pair incompati- 
bilities, each chosen independently at random, and the first and the second member of each pair chosen 
independently from the an available alleles. We will assume that m = \ca^n, label c' = \c, and denote, 
as usual, the resulting set of incompatibilities by F. 

To see that these two models are equivalent for our purposes, first note that the number of incom- 
patibilities which are not legitimate, in the sense that the two alleles are chosen from the same locus, is 
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stochastically bounded in n. (In fact, it converges in distribution to a Poisson(c ^ a^) ra ndom variable.) 
Moreover, by the Poisson approximation to the birthday problem (iBarbour et al. , 1992 ). the number of 



pairs of choices which result in the same incompatibility in this model is asymptotically Poisson(c'a^/2). 
In short, then, the procedure results in the number m — 0(1) of different legitimate incompatibilities. If 
m in the modified model is increased to, say, m' = m + v?!'^ , then the two models could be coupled 
so that the incompatibilities in the original model are included in those in the modified model. As the 
existence of a viable phenotype becomes less likely when m is increased, this demonstrates that (fTTI) will 
follow once we show the following for the modified model: for every e > there exists a large enough 
a so that c' < log a — e implies that > 1 a. a. s. 

To show this, we introduce the auxiliary random variable 

^ = X] n (^ol{|/n<7|=o} +^^il{|/na|=i}) , 

where \a is the indicator of the set A. The size of the intersection / n o" is computed by transforming 
both the incompatibility / and the genotype a to sets of (indexed) alleles, and the weights wq and w\ 
will be chosen later. To intuitively understand the statistic X, note that when w^^ = w\ = 1, the product 
is exactly the indicator of the event that a is viable and X is then the number of viable genotypes A^. In 
general, X gives different scores to different viable genotypes — however, the crucial fact to note is that 
that A > iff A > 0. Therefore 

P(A > 0) = P(X > 0) > {E{X)f lEiX^), 

which is how the second moment method is used ( Achlioptas and MooreL 2004 ). 

As 



p(|(7n/| = 0) 
P(kn/| = 1) 



1 \ 2 

a — 1 



a 

2(a- 1) 



we have 



Moreover 




E^X") = E'^^Q (« - l)'(^o^'(00) + 2woWiP{01) + wlPill)), 

where P(01) is the probabihty that / has intersection of size with a = Oi . . . O^Ofc+i . . . 0„ and of size 
1 with T = li . . . IfcOfc+i . . . On, and P(00) and P{ll) are defined analogously. Thus, if A; = an, 

„/ X / 1 + ^ 
P(00) =1 — 



X 2a / 1 + a 
P(01) = — 1 



a \ a 



P(ll) = ^il^fl-i±i^V2^" 



a \ a / \a 
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Let A = Aa,wo,wi (a) be the n'th root of the k = (an)'th term in the sum for E{X^), divided by E{Xy^ 
Hence 

A- <°-l)° 



a ■ a°(l - 

(1 + ""o".! (1 - + (1 - + (f)^ 

Let a* = (a — l)/a. A short computation shows that A = 1 when a = a*. 

If A > 1 for some a, then E{X'^)/{E(X)y increases exponentially and the method fails (as we 
will see below, this always happens when wo = wi = 1, i .e., when X = N). On the o ther hand, if 
A < 1 for a 7^ a*, and ^(a*) < 0, then Lemma 3 from lAchlioptas and Moord (|2004f) implies that 
E{X^)/{E(X)y < C for"some constant C, which in turn imphes that P{N > 0) > 1/C. The sharp 
threshold result then finishes off the proof of (fTTI) . 

Our aim then is to show that wq and wi can be chosen so that, for c' = log a — e, A has the properties 
described in the above paragraph. We have thus reduced the proof of ([TT]) to a calculus problem. 

Certainly the necessary condition is that ^(a*) = 0, and 

^(«*) = -^{Ma -l)-wi{a- 2))\ 

so we choose wq = a — 2 and wi = a — 1. (Only the quotient between wq and u^i matters, so a single 
equation is enough.) This simplifies A to 



A A / \ [a -I 

A = K{a) 



s2 N 2c' a? 



Let If = log A. We need to demonsti-ate that (/9 < for a G [0, a*) U (a*, 1] and that (p"{a*) < 0. A 
further simplification can be obtained by using x — Cx^ < log(l + x) < x (valid for all nonnegative x), 
which enables us to transform ip (without changing the notation) to 

^4 / ^ 1 \ ^ 

ip(a) = c'- -j- ( a I — Q log a — (1 — q) log(l — a) + a log(a — 1) — log a. 

(a -1)4 V a J 

Now 

1 

ip"{a) = 2c'- 



(a-l)4 a(l-a)' 

So automatically, for c' large but c' = o{a), ip"{a*) < for large a. Moreover, ip cannot have another 
local maximum when (p" > 0. If 99(a) > for some a / a*, then this must happen for an a in one of 
the two intervals [0, l/{2c') + 0{{c'y^)] or [1 - l/(2c') - ©((c')"^), 1]. Now, has a unique maximum 
at a* in the second interval. In the first interval, a short computation shows that 



(p{a) = —e — a log a + O 



log log a 
log a 
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Table 1 : The lower bounds on 7 obtained by the method described in text, compared to the easy upper 
bounds 4 log a. 



a 


1. b. on 7 


4 log a 


3 


1.679 


4.395 


4 


2.841 


5.546 


5 


3.848 


6.438 


6 


4.714 


7.168 


7 


5.467 


7.784 


8 


6.128 


8.318 


9 


6.715 


8.789 


10 


7.242 


9.211 


20 


10.672 


11.983 


30 


12.608 


13.605 


40 


13.944 


14.756 


50 


14.960 


15.649 


100 


18.017 


18.421 


200 


20.982 


21.194 


300 


22.663 


22.816 


400 


23.846 


23.966 


500 


24.759 


24.859 



which is negative for large a. This ends the proof. 

This method yields nontrivial lower bounds for 7 for all a > 3, cf. Table 1. 



Appendix E. Existence of viable phenotypes. 

In this section we describe a comparison between models from Sections 5.2 and 5.3 that will yield the 
result in Section 5.3. We begin by assuming that a = 1/r is an integer, which we can do without loss 
of generality. Divide the i'th coordinate interval [0, 1] into a disjoint intervals Iio, . . . , /i,a-i of length r. 
For a phenotype x G P let A(a;) G Qa be determined so that A{x)i = j iff Xi G lij. 

Note that, as soon as Jj^j^ x li^j^ contains a point in Viji.2, no x with A(.t) jj = ji and A(a;)j2 = J2 is 
viable. This happens independently for each such Cartesian product, with probability 1 — exp(— Ar^) > 
cr^/ (2n). Therefore, using the result from Section 5.2, when cr^ > 4 log a = — 4 log r, there is a. a. s. no 
viable genotype. 

On the other hand, let P be the closed e-neighborhood of the interval 7 in [0, 1] (the set of points 
within e of /), and consider the events that I-^j^ x /^^^ contains a point in Ilj^jj- These events are 
independent if we restrict ji, j2 to even integers. Moreover, each has probability 1 — cxp(— 4Ar^) ~ 
4cr^/(2n), for large n. It again follows from Section 6.2 that a viable genotype x with A{x)i even for 
all i, a. a. s. exists as soon as 4cr^ < 4(log(a/2) — o(l)) = (— 41ogr — log 2 — o(l)). 



